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1.  Introduction 


Tne  report  on  "Applications  of  a  Dedicated  Gravitational 
Satellite  Mission"  (National  Academy  of  Sciences,  1979) 
recommended  rigorous  simulation  analyses  to  answer  questions 
on  the  resolution,  accuracy,  economical  data  reduction 
procedures,  etc. ,  for  an  improved  gravity  field  from  the 
GRAVSAT  mission. 

Several  error  analyses  have  been  reported  recently 
for  both  the  "High-Low"  and  "Low-Low"  missions.  Plsacane 
and  Yionoulis  (1980),  and  Douglas  et  al.  (1980)  used  least 
squares  "deterministic"  techniques  to  get  estimate  of  l*x  1° 
anomaly  recovery  errors.  Lancaster  et  al.  (1980)  used 
similar  techniques  but  obtained  less  optimistic  estimates. 
Least  squares  collocation  techniques  were  applied  by  Rapp 
and  Hajela  (1979)  for  the  "High-Low"  case,  and  by  Rummel 
(1980)  for  the  "Low-Low"  case.  Krynski  (1979)  e^so  used 
the  collocation  method  for  his  error  analysis.  Jekeli 
and  Rapp  (1980),  Breakwell  (1980),  Breakwell  and  Duhamel 
(1981)  related  the  spectrum  of  the  measurement  of  range- 
rate  and  gradiometry  to  the  spectrum  of  the  disturbing 
potential  of  the  earth’s  gravity  field  to  get  better  insight 
into  resolution  and  achievable  accuracies. 

These  error  analyses  have  helped  to  define  the  current 
GRAVSAT  mission  configuration  plan  (Bergeson-Villis  et 
al.,  1981)  of  two  low  satellites  at  an  altitude  of  160 
km,  in  the  same  polar  orbit  but  separated  by  300  km.  The 
uncertainty  of  relative  range-rate  measurement  is  expected 
to  be  less  than  1  ym/s.  Wong  and  Sjogren  (1980)  describe 
the  simulation  of  range-rate  data  from  1500  30’x30'mean 
anomalies  residual  to  GEM9  (Lerch  et  al.,  1979)  for  a  config¬ 
uration  closely  corresponding  to  the  current  GRAVSAT  config¬ 
uration.  The  anomalies  are  in  and  around  the  sharply  varying 
field  associated  with  the  Tonga  Kermedec  Trench  in  the 
south-west  Pacific  Ocean.  This  data  is  described  in  Section 
3  and  will  be  used  in  the  present  simulation  study. 

The  least  squares  collocation  method  is  used,  described 
in  Section  2,  to  predict  the  30 'x  30' ,  l°x  1®  and  2°x  2*  residual 
anomalies.  The  predicted  anomalies  are  compared  with  the 
input  anomalies  used  for  generating  the  data.  The  purpose 
of  the  present  study  is  thus  to  test  the  prediction  of 
gravity  anomalies  from  the  GRAVSAT  mission,  to  draw  conclu¬ 
sions  for  the  data  requirement  for  economical  (in  terms 
of  computer  time)  prediction,  and  to  establish  the  resolution 
achievable  in  terms  of  the  block  size  of  the  mean  anomalies. 

The  computation  of  line  of  sight  accelerations  sensed 
at  the  satellite  altitude  is  described  in  Section  4.  The 
computation  of  the  covariances  needed  in  the  collocation 
method  is  discussed  in  Section  5  while  the  tests  for  the 
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effect  of  the  extent  and  density  of  data  on  the  prodie tod 
anomaly  aro  diacuoeod  in  Section  6.  The  numerical  rooults 
with  the  resulting  modest  data  requirement  in  terms  of 
computer  time  are  discussed  in  Section  7  for  the  prediction 
of  1°  anomalies.  The  resolution  achievable  in  terms  of  anom¬ 
aly  block  sise  is  discussed  in  Sections  S  and  9  by  considering 
the  prediction  of  30*  and  2*  anomalies  respectively.  A 
representative  sample  of  anomalies  is  chosen  in  a  compar¬ 
atively  smooth  area  away  from  but  near  the  trench ,  and  the 
sample  is  then  extended  to  include  anomalies  in  the  area 
approaching  and  then  directly  over  the  trench. 


2.  Mathematical  Model 


We  define  the  residual  range-rate  p  between  the 
two  satellites  in  the  earth's  anomalous  potential  field  T 

(2.1)  p  -  p0-  pc 


(2.2)  T  -  V  -  U 


where  p0  is  the  range-raje  in  the  earth's  gravi tat ions  1 
potential  field  V  and  pc  is  the  range-rate  in  the 
earth's  normal  potential  field  U  ,  which  was  taken  in 
this  study  to  be  described  by  GEM9. 

It  is  easily  seen  that 

(2.3)  p#  *  kxl  •  eia  ,  xla  ■  Xj  -  Xj 


where  x  ,  x  are  the  position  and  velocity  vectors  of 
the  satellites  in  a  geocentric  quasi- inertial  cartesian 
coordinate  system;  subscripts  1  and  2  refer  to  the  lead 
and  trailing  satellites  respectively  e.g.  in  a  north  to  south 
arc  or  vice  versa;  and  eu  is  the  unit  line  of  sight 
vector  from  the  lead  to  tEe  trailing  satellite.  We  assume 
that  the  range-rate  is  measured  at  the  lead  satellite. 
Equation  (2.3)  expresses  the  fact  that  the  range-rate  is 
the  projection  of  the  relative  velocity  of  the  satellites 
on  the  line  of  sight  between  the  satellites. 

The  time-derivative  of  the  range-rate,  which  we  may 
call  as  line  of  sight  acceleration,  is  obtained  from  (2.3) 
as: 

(2.4)  p0  «  *  e12  +  xXI  •  eia 


and  we  define  the  residual  line  of  sight  acceleration  £ 
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analogously  to  (2.1)  as: 


< 2.5 )  p  -  (Jo  -  pc 

Rummel  (1980,  p,5,  Table  1)  hae  shown  that  the  second 
term  in  equation  (2.4)  is  negligibly  small  as  compared  to 
the  first  term  when  applied  for  the  residual  line  of 
sight  acceleration  p  ,  and  we  may  write  with  a  negligibly 
small  error  (ibid.,  p.ll) 

(2.6)  p  a  k12  *ei2»(x2  -  xx)  •  ei2  *  <22.2  -  23U>  * 


where  VT  is  the  gradient  of  the  anomalous  potential  at 
the  satellite  location.  Again,  the  assumption  in  (2.6) 
is  that  p  is  sensed  at  the  location  of  lead  satellite 
analogously  to  the  measurement  of  range-rate  at  that  location 
in  (2.3). 

We  will  henceforward  treat  p  as  the  "observation" 
related  to  the  anomalous  potential  at  satellite  location 
through  (2.6).  The  reduction  01  p  from  p  will  be  des¬ 
cribed  in  Section  4. 

We  may  denote  VT  in  terms  of  its  components  <$r  ,  <5^,  * 

6^  in  a  geocentric  earth- fixed  spherical  coordinate  \.r,<J>'>*) 
system.  Then  following  Heiskanen  and  Moritz  (1967,  Section  6-4), 
we  have: 

(2.7)  7T  -  8r  er+St,  +  «x  ex 

where  6  „  85  3T/3r 
r 

,=  (1/r)  (3T/3V)  8  *y  5 

<P 

5  *  (1/r  cos  <t>' )  (3T/3X)  *  -yn 

A 

where  y  is  the  normal  gravity  at  the  satellite  location 
and  5  ,  n  are  the  latitudinal  and  meridional  components 
of  the  deflection  of  vertical  at  the  satellite  location. 

We  use  ax ,  bx,  clf  to  denote  the  direction  cosines 
of  line  of  sight  (lead  satellite  1  to  trailing  satellite  2) 
wich  r  ,  <j)'  ,  X  axes  at  location  of  lead  satellite; 
and  similarly  a2  ,  b2  ,  c2  as  direction  cosines  of 
line  v.f  sight  (trailing  satellite  2  to  lead  satellite  1) 
with  1  ,  <j>'  ,  X  axes  at  location  of  trailing  satellite. 

Then  by  inserting  (2.7)  into  (2.6),  we  get: 

(2.8)  p  *  ( VT2  -  VT  1 )  •  e12  **  -VV2  •  621  -  VTi  *  £12 

-  -  a 26r2  +  b2v2£2  +  c2y2r2  -  axfiri  +  biYiSi  +  CiYihi 
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The  direction  cosines  a  ,  b  ,  c  of  the  line  of  sight 
may  be  rigorously  computed  by  relating  the  geocentric 
cartesian  and  spherical  coordinates  cf  the  satellites,  e.g. 
see  Hajela  (1979,  equations  (10)  to  (13)).  For  the  simple  config 
uration  of  satellites  in  this  simulation  study,  the  direction 
cosir.es  may  also  be  obtained  from  simple  geometric  consid¬ 
erations  (see  Section  3). 

The  predicted  value  Ag'  of  a  mean  residual  gravity 
anomaly  is  obtained  from  p  5  T^  (for  ease  of  notation) 
"observations"  by  least  squares  collocation  (Moritz,  1980, 
Section  14): 

(2.9)  Ag'  -  T£  ,  C*  -  £?TJl,TJt  +  - 


-AT 


■A 


where  C^^TS,  -  C0V^Zi,  »  T^)  is  the  covariance  matrix  of 
residual  accelerations  £  =  T&  ;  D  -  a.f.^1,  is  noise  matrix 
of  Tjj,  ,  taken  as  a  diagonal  matrix'Vith  equal  variance  (cr Tg,) 
of  Tj,  data  on  the  diagonal  (see  Section  4);  £Ag'TJl  2  cov(Ag'  ,T^) 
is  the  covariance  vector  of  7$,  with  the  predicted  mean  anomaly 
block,  and  the  superscript  T~"  is  for  the  transpose.  The 
covariance  with  the  mean  block  is  obtained  by  numerical  inte¬ 
gration  of  the  point  covariances  after  dividing  the  block 
into  mxm  sub-blocks.  By  numerical  tests,  it  was  found 
that  m*4,  5,6  is  adequate  for  30' ,  1°  and  2°  blocks.  The 
estimating  operator  A  and  the  underlying  covariances  cov(T^,T£) 
and  cov(Ag',T£)  will  be  discussed  in  Section  5. 

The  residual  anomaly  Ag'  is  obtained  from  anomaly  Ag  , 
referred  to  an  ellipsoidal  field,  by  subtracting  the  anomaly 
Agr,  implied  by  the  normal  potential  U  ,  e.g.  see  Rummel 
et  al.  (1976,  p.20). 


The  estimate  of  the  variance  a2^g'  of  the  predicted 
residual  anomaly  is  given  by: 


(2.10)  32gt 


C0  -  C 


c*-i  c 

Ag'TA  -  ^Ag'TA 


where  C0  is  the  variance  of  the  residual  anomalies  of  specified 
block  size  and  is  dependent  on  its  latitude  in  the  case  of 
equi-angular  blocks. 

The  covariances  needed  in  (2.9)  and  (2.10)  are  computed 
by  the  propogation  of  covariances  utilizing  (2.8).  For  example, 

for  the  data  points  i  and  j;  i  -  1  ,  . . . ,  n  ,  j»l .  i  i 

the  corresponding  element  in  the  matrix  C?*  is  obtained 
by:  *  * 


The  computation  of  (2.11)  and  (2.12)  is  discussed  further 
in  Section  5.  Some  simplification  also  results  due  to  similar 
values  of  some  direction  cosines  in  this  simulation  study. 

The  auto-  and  cross-covariances  of  <Sr  ,  £  ,  n  ,  Ag'  are  computed 
from  subroutine  COVAX  (Tscherning,  1976)  after  modifying 
it  to  give  all  required  covariances  in  mgals2  to  correspond 
to  the  units  of  (Tg^T^j)  and  (Tgi,Ag')  .  The  degree  variances 
for  Z «  3  , . . . ,  20  were  made  zero  as  we  are  considering  residual 
accelerations  p  =  Tg  and  anomalies  residual  to  G1M9  ,  which 
is  complete  to  degree  and  order  20.  Though  the  low  degree  var¬ 
iances  in  GEM 9  and  the  model  are  strictly  not  the  same, 
this  approximation  has  negligibly  small  effect  on  equations 
(2.9)  and  (2.10).  The  underlying  anomaly  degree  variance 
(c^)  model  used  is  (see  Tscherning  and  Rapp,  1974): 

(2.13)  c^  -  mgals2  ;  r0  -  6,369,779.8  m  , 
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where  rg  is  the  radius  of  the  Bjerhammar  sphere  internal 
to  the  earth's  surface,  so  that  the  degree  variance  cr^CAg) 
at  a  point  with  geocentric  distance  r  is  given  by : 

(2.14)  o^(Ag)  ■  s*+2  cA  ;  s  ■  (rg/r)2  . 


3.  Residual  Range  Rate  Data 


Wong  and  Sjogren  (1980,  Sec. 2. 2)  simulate!  the  range 
rate  between  two  low  satellites  at  an  altitude  of  150  km  above 
a  spherical  earth.  The  two  satellites  were  in  the  same  polar 
orbit,  i.e.  same  geocentric  longitude  for  both  satellites 
for  any  range  rate  measurement,  with  a  separation  of  200  km  * 
1“7554. ..  between  the  two  satellites.  The  satellite  ephemeris 
was  perturbed  by  the  anomalous  gravity  field  defined  by  1500 
30'  x30'  anomalies  covering  the  area  bounded  by  latitudes 
0°  to  25°  S.  and  longitudes  180°  to  195°  E.  These  anomalies 
were  obtained  from  altimeter  measurements  (R.H.  Rapp,  private 
communication,  1980),  and  from  these  the  anomaly  implied  by 
GEM9  was  subtracted  out.  Accordingly,  the # simulated  measure¬ 
ment  corresponded  to  residual  range  rate  p  in  (2.1). 

A  random  noise  (mean  *  0  ,  o  *»  lym/s)  was  impressed  on 
p  to  simulate  the  precision  of  the  GRAVSAT  mission.  The 
integration  time  and  the  data  interval  along  an  arc  was  2 
seconds  ( A4>  s  0?14),  while  the  separation  between  arcs  in 
longitude,  AX  ,  was  0?2  in  the  central  half  of  the  area 
and  0.4°  in  longitude  on  each  side.  This  density  of  arcs 
was  obtained  for  both  N/S  (lead  satellite  moving  from  north 
to  south)  and  S/N  (lead  satellite  moving  from  south  to  north) 
arcs.  A  total  of  23,128  "noisy"  £  measurements  were  generated 
over  59  each  N/S  and  S/N  arcs,  and  sent  to  us  on  a  magnetic 
tape  for  this  study  (L.  Wong,  private  communication,  1980). 

A  portion  of  the  input  30'  residual  anomalies  are  shown 
in  Figure  3.1.  These  have  been  contoured  at  20  mgals  interval 
in  Figure  3.2  to  show  the  gradient  of  the  anomaly  field.  We 
notice  that  the  field  is  comparatively  smooth  from  the  north 
edge  up  to  12°  S.,  varies  sharply  as  the  trench  is  approached 
around  14°  S.  from  where  the  trench  proceeds  south  generally 
about  the  longitude  187?5  E  .  The  area  for  testing  the  prediction 
of  l°x  1°  anomalies  was  accordingly  selected  between  8®  to  16°  S. 
and  186°  to  188°  E.  (and  slightly  to  east  and  west)  astride 
the  steep  gradients  of  the  trench,  and  also  covering  a  compar¬ 
atively  smoother  field  in  the  northern  portion .  The  1°  anom¬ 
alies,  (which were  predicted  based  on  $  generated  from  the 
input  30'  anomalies,)  reflect  the  general  gradient  of  the  30' 
anomalies.  The  1°  anomaly  f^eld  will  be  discussed  in  Sections 
6  and  7. 

The  p  measurements  on  the  data  tape  were  identified 
with  time  along  the  arc  from  specified  starting  position  of 
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Figure  3.1  A  Portion  of  the  Input  30'x30'Mean  Residual 
(to  GEM9)  Anomalies  in  mgals.  1500  anomalies  covering 
0°  £  <J>  <  25°  S.  ,  180°  £  X  <  195°  E.  were  used  to  simulate  0  data 
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ttr  lead  (and  trailing  satellite)  but  the  actual  location  of 
the  satellies  for  a  particular  P  was  not  given,  except  for 
sample  runs  for  one  N/S  and  S/N  arc  at  intervals  of  3C  seconds. 
The  geocentric  spherical  coordinates  (r,<J>',X)  of  the  satellites 
were  computed  by  us  in  the  following  simple  manner  and  agreed 
with  the  sample  values  to  at  put  0V004  and  O'.'OOl  in  latitude 
and  longitude  respectively.  As  the  satellite  orbits  were 
generated  over  a  spherical  earth,  the  period  P  of  the  satellite 
is  given  by 

(3.1)  P  *  2ir  a3/i  /(GM)*  ,  a  -  ae  +  150  km 

a  5249  s  for  the  specified  values  of  earth  radius 
ae  and  the  gravitational  constant  times  earth's  mass  GM  . 

Further,  as  the  satellties  are  in  a  circular  orbit  and  the 
inertial  starting  azimuth  is  180°  or  0°  for  N/S  and  S/N  arcs 
respectively,  the  change  in  declination  (A6)  or  geocentric 
lat  ude  (A<J>')  was  constant  over  short  arc  span  of  about 
6i  luinutes;  and  we  have  for  the  data  interval  of  2  seconds? 

OCA  O 

(3.2)  A$'a  A5  »  •  2  with  A<{> '  in  degrees  and  P  in 

p  seconds. 


The  geocentric  longitude  of  the  satellite  changes  only 
due  to  the  rotation  of  the  earth,  and  finally  the  height  h 
above  the  ellipsoid  may  be  calculated  from: 

(3.3)  h»a-p  ,  p  -  a(l  -  e2)*/(l  -  e2  cos2  <*>'  )* 


Here  p  is  the  radius  vector  to  the  ellipsoid  at  geocentric 
latitude  $ '  ,  and  e  is  the  first  eccentricity.  However, 
as  the  satellite  orbits  were  generated  over  a  spherical  earth, 
h  was  uniformly  taken  as  150  km. 

The  direction  cosines  in  equation  (2.8/  of  the  line  of 
sight  with  the  (r,<J>',A)  axes  at  the  two  satellites  (k*l,2)  for 
different  p  (and  p)  data  points  ( i *  1 , . . .  ,  n)  can  be  simply 
obtained  from  geometric  considerations  because  of  the  simple 
configuration  of  satellites.  This  is  shown  in  Figures  3.3a 
and  3.3b,  where  two  data  points  i  ,  i+1  have  been  considered 

for  a  N/S  and  a  S/N  arc  respectively.  As  the  two  satellites 
are  in  a  circular  orbit  the  line  of  sight  %k  makes  equal 
angles  with  the  radial  direction  at  k»l,2  for  all  i  , 
and  this  is  a  constant  angle  because  of  constant  separation 
between  the  satellites  of  200  km  ( 9  a  1?7554 . . . a  1°  45*  19"). 

Also  the  longitude  of  the  two  satellites  is  the  same  for  all 
i  .  Accordingly,  the  angles  aik  ,  0^  ,  Yik  of  £ik  with 
the  positive  direction  of  axes,  and  thus  the  direction  cosines 
aik  »  bik  *  °ik  are  as  follows: 
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The  simulated  residual  range  rate  ft  was  first  filtered  j 

for  the  "observational"  noise  of  lym/s  by  approximating  the 
ft  values  in  the  least  squares  sense  with  a  cubic  spline  function 
(fts)with  fixed  knots  (see  Hajela,  1979,  equations  (17)  to 
(19)).  The  spline  function  (5  was  then  differentiated  to  • 

get  p  .  The  degree  of  smoothing  or  the  filtering  of  noise 
depends  on  the  spacing  of  the  spline  knots;  a  larger  spacing  j 

results  in  a  greater  smoothing.  Several  different  knot  spacing  j 

were  tested  at  4,  6,  8,  10,  12,  14,  IS  seconds  for  a  maximum 
of  121  p  data  points  at  2  seconds  interval  in  several  N/S 
and  s/N  arcs.  The  difference,  Aft  ,  between  the  smoothed  f 

spline  function  fts  and  the  noisy  ft  data  was  examined  after  1 

discarding  the  data  points  in  two  intervals  each  in  the  beginning  j 

and  the  end  of  the  approximating  spline  function.  This  was 
done  to  discard  the  data  affected  by  spurious  oscillations  1 

in  the  spline  because  of  arbitrary  end  conditions. 

(4.1)  Aft  -  ftg  -  ft  ; 


Some  statistics  for  one  N/S  and  one  S/N  arc  in  the  center 
of  the  area  are  given  in  Table  4.1. 


Table  4.1  Filtering  of  Noise  in  P  for 
Different  Spline  Knot  Spacing 


SALINE  KNOT 

SPACING 
( sec . ) 

# 

PTS . 

RMS  VALUE 

&  a 

(ym/s)  8 

STATISTICS  FOR  Ap 

MEAN  RMS  Min . 

(Um/8) 

Max. 

N/S  ARC 

4 

47 

164.5 

164.5 

0.01 

0.92 

-2.28 

1.35 

6 

70 

164.9 

164.8 

-0.01 

0.91 

-2.17 

1.50 

8 

93 

234.5 

234.5 

-0.02 

0.92 

-2.46 

1.52 

10 

101 

256.1 

256.1 

-0.01 

0.98 

-2.32 

1.77 

12 

97 

244.8 

244.8 

-0.04 

1.41 

-2.68 

3.10 

14 

81 

229.5 

229.4 

0.00 

2.32 

-6.96 

6.26 

16 

89 

224.3 

224.3 

0.14 

5.29 

-13.30 

14.75 

S/N  ARC 

4 

47 

133.6 

133.6 

-0.01 

0.89 

-2.28 

1,36 

6 

70 

169.0 

169.0 

0.00 

0.93 

-2.43 

1.59 

8 

93 

163.1 

163.1 

0.00 

0.90 

-2.54 

1.56 

10 

101 

180.6 

180.6 

-0.01 

0.^7 

-2.61 

1.98 

12 

97 

170.6 

170.5 

0.00 

1.17 

-3.42 

3.14 

14 

91 

160.4 

160.3 

-0.03 

1.99 

-5.04 

4.94 

16 

89 

158.6 

158.4 

-0.03 

3.99 

-8.93 

10.78 

We  note  from  Table  4.1  that  to  filter  the  noise  of  mean  0 
and  standard  deviation  lym/s,  spline  huot  spacing  of  10 
seconds  appears  to  be  suitable.  The  spline  function  with 
this  spacing  was  accordingly  used  co  approximate  p  in  all 
arcs.  For  101  data  points  in  each  arc,  mean  of  ranged 

from  >0.01  to  0.00um/s  and  RMS  &p  ranged  from  0.92  to 
1.04ym/s  .  The  first  derivative  of  the  spline  function  gave 
the  p  values.  The  $  and  p  values  for  one  S/N  arc  over 
the  nor..-'  '.al  longitude  of  187°E.  ar«  shown  in  Figures  4.1  and 
<i.2  respectively.  The  variation  in  the  H>  and  p  values 
will  be  discussed  in  Section  5. 

Ruromel  (1980,pp.  13,  43)  has  discussed  the  corresponding 
precision  of  p  based  on  the  precision  of  p  by  comparison 
of  the  degree  at  which  the  signal  to  noise  ratio  becomes  one 
in  the  two  cases.  For  the  present  study  with  satellite  alti¬ 
tude  of  150  km,  the  above  reasoning  gives: 

(4.2)  0^  =  oT4  ■  0.05  mgals  for  a ^  ■  lum/s  at  h  -  150  km 


Accordingly,  0.05  mgals  will  be  used  for  oTo  in  equation 
(2.9),  where  D  -  <j*Tz  I  .  * 
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Figure  4.1  Smoothed  Residual  Range  Rate  (wn/s)  in  a  Central  S/M  Arc 

♦'  -2151  to  -456  h  »  150  km 
X  18753  to  18653  9  «  200  km 


Residual 


S»  Auto-  and  Cross-Covariance  Functions 


We  need  the  covariances  cov(T^,T4)  and  cov(Ag'  ,T.) 
for  predicting  a  residual  anomaly  Ag”  in  equation  (2.3): 

<2.9)  Lf  ■  I,  .  i*  -  £«,«♦  fi  •  fi  *  ®T»  i 

-ill 


where  A  is  the  estimating  operator  indicating  the  sensitivity 
of  X&  "data  in  predicting  the  anomaly,  or  a  weighting  function 
(Breakwell,  1981),  which  may  be  examined  to  judge  if  distant 
data  may  be  excluded  without  significantly  affecting  the^predic- 
ted  value. 

We  will  accordingly  examine  the  underlying  oovariances 
and  the  operator  A  to  see  if  a  modest  extent  and  density 
of  data,  which  will  be  economical  in  computer  time,  will  lead 
to  a  "satisfactory"  prediction.  The  latter  may  be  checked 
by  the  discrepancy,  e  ,  of  prediction  A* g'  compared  to  the 
input  anomaly  Ag'  used  or  implied  in  simulating  the  data: 

(5.1)  e  »  Ag'  -  Ag' 


The  computation  of  the  covariances  in  (2.11)  and  (2.12) 
is  somewhat  simplified  for  the  present  study  in  view  of  Figure 
3.3  and  equation  (3.4)  for  the  direction  cosines.  We  recall 
that  for  all  data  points  i  in  N/S  arcs,  or  the  S/N  arcs, 
the  direction  cosines  at  the  two  satellites  1,2  are: 


ail 

"  ai2 

an 

-c 

(5.2)  bu 

"  “bi2 

m 

+0 

bil 

*  -bi2 

- 

-0 

cil 

"  ci2 

m 

0 

015,318. . .  •  -a 

for 

999,882...  •  b 

for 

999,882. . .  »  -b 

for 

for 

both  N/S  and  S/N  arcs 
N/S  arcs 
S/N  arcs 

both  N/S  and  S/N  arcs 


Accordingly  the  first  vector  transposed  in  (2.11)  and 

(2.12)  becomes:  [a  ,  -byjo  »  a  »  bYn]  for  N/S  arcs,  and 

the  sign  of  b  is  reversed  for  S/N  arcs.  Also,  while  modifying 
the  units  of  covariances  to  be  in  mgals*  for  COVAX  (see  below 

(2.12) )  if  we  multiply  the  covariances  involving  €  and  n  by 
-y  ,  then  using  (2.7);  equation  (2.11)  is  changed  to: 


6rj2  VJ2  Vj1 


(5.3) 

<firi2,6rj2)  *  *  ^riR’^jl* 

r  - 

a 

<fi*»i2»6rj2)  *  *  (fi^li2»6Wl) 

b 

e”<TU*V  ' 

<6ril,6rj2)  ‘  *  (,Sril*6^jl) 

a 

«  <TU.T„> 

*  *  (S4'il»6^jl> 

-b 

«  [a  ,  b  ,  a  ,  -b]  • 


Similarly  equation  (2.12)  is  changed  to: 


(5.4) 

5ri2*A«‘ 

cov(Tai,Ag') 

^•ia**8’ 

6ril»*8' 

*  [a  ,  b  ,  a  ,  -b]  • 

mm*  mmm 

In  (5.3)  and  (5.4',  the  sig-  of  b  applies  for  N/3  arcs; 
and  it  should  be  reversed  for  S/N  arcs. 

With  two  satellites  each  at  points  i  and  j  ,  there 
are  six  distances  (ty)  to  be  considered  for  the  computation 
of  covariances  in  (5,3).  If  j  ■  i  ,  the  distance  is  either 
0  or  0  ,  the  separation  between  the  satellites.  However, 
when  j  #  i  ,  the  distances  il  -  jl  ,  i2  -  j2  ,  il  -  j2  , 
i2  -  jl  are,  in  general,  all  different. 

For  n  data  points,  n(n+l)/2  covariances  are  required 
in  the  matrix  and  each  element  requires  16  calls  to 

COVAX  for  the  commutation  of  (5.3).  For  the  tests  in  this 
study,  the  covariances  were  rigorously  computed  to  avoid  any 
approximations  in  this  regard.  However,  for  production  runs 
in  future,  interpolation  from  covariance  elements  stored  in 
a  table  would  be  faster,  e.g.  see  Sttnkel(19?9) . 

The  cov<Tm,T£j)  function  has  been  plotted  in  Figure 
5.1.  This  refers  to  residual  line  of  sight  accelerations 
at  150  km  altitude,  residual  to  degree  20,  and  a  separation 
of  200  km  between  the  satellites.  The  graph  is  comparable 
to  the  second  curve  in  Figure  5a  cf  Rummel  (1980,  p.19). 
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Figure  5.1  Covariance  Function  of  Residual  (to  GEM9) 

Line  of  Sight  Acceleration.  Satellite  Altitude 
150  km,  Satellite  Separation  200  km 
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The  correlation  length  in  our  case  is  1?3  so  that  the  predicted 
anomaly  would  be  most  sensitive  to  data  within  \p  ■  1?3 

(Moritz.  I960,  p.177).  The  extent  of  data  will  be  numerically 
tested  in  Section  6. 

The  cov(Tni,Ag')  function  for  one  of  the  mean  l°x  1° 
residual  anomalies  has  been  plotted  in  Figure  5.2.  Similar 
graphs  are  obtained  for  other  anomalies  also.  The  two  graphs 
in  Figure  5.2  show  the  covariance  of  anomaly  block  -13P>  <f>  >  -14° 
187°  <  X  £  188°  with  a  N/S  arc  and  a  S/N  arc  each  with  nominal 
longitude  187?5. 

We  note  that  the  correlation  length  is  1?2;  and  that 
the  peak  of  the  covariance  function  does  not  occur  when  the 
lead  satellite  is  over  the  center  of  the  anomaly  block*  but 
when  it  is  past  the  center  by  a  distance  of  0/2  ,  where  e  is 
the  separation  between  the  satellites. 

The  location  of  the  lead  and  trailing  satellites  and 
the  anomaly  block  are  specifically  accounted  for  in  (5.3) 
and  (5.4),  and  the  location  of  the  measurement  is  specifically 
accounted  for  in  (2.3)  and  (2.3).  However,  the  graphs  in 
Figure  5.2  imply  that  the  selection  of  T^  data  (sensed  at 
the  lead  satellite)  should  be  centered  (e.g.  for  a  spherical 
distance  ^  >-  1°  sc  correlation  length)  around  a  point  which 
is  down  the  arc  at  a  distance  0/2  from  the  anomaly  block 
center.  This  agrees  with  Breakwell's  interpretation  (1981, 
p.10)  of  a  data  "point"  as  the  mid-point  between  the  two  satl- 
lites  at  the  time  of  a  measurement. 

*.Y a  may  also  interpret  this  intuitively  as  shown  in  Figure 
5.3.  Unless  the  separation  between  the  satellites  is  "small" 
compared  to  the  anomaly  block  size,  the  lead  satellite,  or 
the  trailing  satellite,  is  primarily  perturbed  (Fig.  5.3a) 
when  it  is  directly  over  the  anomaly  block.  The  relative 
range  rate  between  the  satellites  varies  at  these  two  instances 
and  hence  a  significant  line  of  sight  acceleration  is  sensed 
(Fig.  5.3b)  at  the  lead  satellite.  Unless  the  separation 
between  the  two  satellites  is  "large"  compared  to  the  anomaly, 
block  size,  both  satellies  are  perturbed  between  the  two  extreme 
locations  in  Figure  5.3a,  and  accordingly  the  data  selection 
should  be  centered  around  a  point  midway  between  these  locations 
as  shown  in  Figure  5.3b.  We  have  implicitly  assumed  here 
that  there  is  an  optimum  separation  distance  compared  to  the 
anomaly  block  size  of  interest. 

If  we  now  wish  to  compare  the  data  points  in  a  N/S  and 
S/N  arc,  both  passing  over  the  center  of  the  anomaly  block, 
to  examine  which  data  points  are  primarily  affected  by  the 
anomaly  block,  then  we  may  conceptually  consider  the  "data 
point"  to  be  midway  between  the  two  satellites.  This  is  shown 
in  Figure  5.3c. 


- - ti- - r-"“? 
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Figure  5.2  Cross-Covariance  Function  Between  Surface  Gravity 
Anomaly  &  Line  of  Sight  Acceleration,  Satellite 
Altitfude  150  km.  Satellite  Separation  200km  (-e). 
Anomaly  and  Acceleration  Residual  to  GEM9 
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S/N  ARC 


Fig.  5.3a 


Nominal  Locations  -  Only  1  Satellite  Primarily 
Perturbed  by  the  Anomaly  Block 


Fig.  5.3b 


Data  Selection  for  Predicting 
Anomaly  Block 

Location  of  Lead  Satellite 
as  in  Fig.  5.3a 


5* 


12 


Figure  5.3  Data  Selection  fot  Predicting  Anomaly  Block  & 
Comparison  of  "Data  Points"  for  N/S  k  S/N  Arcs 
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We  should  again  point  out  that  the  location  of  data 
points  remains  at  the  location  of  lead  satellite  in  accordance 
with  (2.6)  and  (2.3)  and  the  actual  location  of  the  two 
satellites  is  used  to  compute  the  covariances  in  (5.3) 
and  (5.4).  However,  the  selection  of  data  points  for  predicting 
an  anomaly  block  is  done  for  a  specified  i|/  ,  as  shown  in 
Figure  5.2.  It  is  only  for  comparison  of  data  points  in 
N/S  arcs,  as  compared  to  S/N  arcs,  to  see  the  effect  of 
anomaly  blocks  on  these  data  points  that  we  need  to  consider 
the  "data  point"  as  midway  between  the  two  satellites . 

The  comparison  of  p  and  p  data  over  N/S  and  S/N 
arcs  is  shown  in  Figures  5.4a,b  and  5.5a,b  respectively. 

First  of  all,  we  notice  from  Figure  5.5  that  t;;ere  is  no 
discernible  difference  at  few  tenths  of  a  mgal  level  between 
the  accelerations  caused  by  the  anomaly  blocks  whether 
it  is  for  N/S  or  S/N  arcs  .  This  is  as  it  should  be.  "fhe 
variation  of  acceleration  due  to  the  gradient  of  the  anomaly 
field  will  be  considered  in  Section  7. 

However,  we  note  from  Figure  5.4  that  though  the  broad 
pattern  of  range  rate  (which  in  fact  reflects  the  variation 
of  p  ,  i.e.,  p)  is  similar,  there  are  very  large  difference 
of  range  rates  in  N/S  and  S/N  arcs.  It  varies  over  the 
area  from  -1.5  to  +1.9mm(s.  It  is  not  clear  if  this  variation 
is  only  due  to  the  satellites  in  the  N/S  and  S/N  arcs  traversing 
different  gravity  field  before  they  reach  a  common  location 
(for  comparison),  and  we  should  expect  this  variation. 

Or  is  it  also  due  to  orbital  errors  introduced  in  the  arcs. 

Wong  and  Sjogren  ( 1980, p. 14)  describe  that  perturbations 
of  the  order  of  0.5  m  were  added  to  some  arcs  in  the  radial 
direction  which  could  cause  perturbation  in  p  of  the 
order  of  -630  to  -930pm/s  over  a  span  of  25°.  As  this 
span  is  traversed  in  6$min. (“390  s  ) ,  a  variation  of  300ym/s 
in  p  will  result  in  a  variation  of  only  0.08  mgals  in  p 
This  only  reflects  the  fact  that  orbital  errors  cause  long 
wavelength  perturbations  in  (5  and  these  are  largely  filtered 
out  in  computing  p  .  (For  some  numerical  results,  see 
Hajela,  1977,  p.64).  It  could  not  be  confirmed  if  the 
data  tape  received  by  us  had  any  orbital  errors  introduced 
or  not. 


We  finally  consider  the  estimating  operator  or  weighting 
function  A  -  CAe,TJl{'*”1  in  e<luation  (2.9).  Figure  5.6a,b 
show  a  plot  of  .Sis^function  for  a  N/S  and  a  S/N  arc  passing 
centrally  over  the  l°x  1°  anomaly  block  -13°  i  <j>  i  -14°, 

187°  <  X  <  188°  (cf  Fig.  5.2).  We  notice  again  that  the 
peak  of  the  function  is  displaced  down  the  arc  from  the 
center  of  anomaly  block  by  one-balf  the  separation  between 
the  satellites.  Secondly,  the  prediction  of  anomaly  block 
is  most  sensitive  to  data  points  within  ^  —  1°  from  this 
peak.  This  interval  has  been  marked  on  the  zero  lines 
of  Figure  5.6a,b.  Thirdly,  the  weighting  function  fluctuates 
with  decreasing  amplitude  as  we  go  farther  from  the  center. 
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Figure  5.4a  Residual  (to  GEM9)  Range  Rate  Data  in  N/S  Arcs 
Contour  Interval  100um/s 

Satellite  Altitude  150  km,  Satellite  Separation  200  kzn 
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Fig.  5.6b  S/N  Are 


Figure  5.6  Estimating  Operator  or  Weighting  Function  for  l°x  1° 
Mean  Residual  Anomalies  from  Residual  (to  GEM9 ) 

Line  of  Sight  Accelerations  (o  m«0.05  mgals) 

Satellite  Altitude  150  km ,  Satellite  Separation  200  km 


Note:  The  dashed  line  is  at  0/2  from  center  of  anomaly  block 
(see  p.21  last  para.) 


This' ref lects  not  only  the  instability  of  the  weighting 
function  contributed  by  C  “l  ,  but  also  the  fact  that  besides 
the  contribution  of  farther  points  being  smaller,  it  also 

+?  largely  cancel  out.  It  should  therefore  be  adequate 
for  "satisfactory"  prediction,  see  equation  (5.1),  to  restrict 
the  extent  of  data  to  ^  ■  1®  .  This  will  lead  to  a  modest 
requirement  of  data  resulting  in  economy  in  computer  time. 

The  instability,  evidenced  by  large  fluctuations  of 
the  weighting  function  was  examined  for  different  noise 
levels,  i.e.,  (see  equation  4.2).  The  weighting  function 

had  increasingly  larger  fluctuations  as  Om^  was  decreased 
to  0.03;  0.02  and  0.01  mgals.  This  is  obviously  due  to 
poor  inversion  of  C*“’  as  with  the  decrease  of  diagonal 
elements  of  D  in  (2.9),  C*  does  not  remain  diagonally 
dominant.  On  the  other  hand,  with  increase  of  to 

0.10,  0.15,  0.20  mgals  the  weighting  function  gets  increasingly 
damped.  We  have  thus  a  numerical  confirmation  of  (4.2), 
i.e.  the  standard  deviation  of  T&  being  0.05  mgals  for 
the  data  used  by  us.  We  have  to  also  ensure  the  symmetry 
of  data  and  a  balance  between  the  density  of  data  and 
the  number  of  points  used  for  the  stability  of  the  estimating 
operator,  e.g.  see  Rummel  et  al.  (1979,  p.355).  Finally, 
it  is  interesting  to  compare  Figure  5.6  with  the  plot 
of  theoretical  weighting  function  in  Breakwell  (1981), 
which  shows  similar  behavior. 

6.  Numerical  Tests  for  Data  Configuration 


The  30'x30'mean  residual  anomalies  in  Figure  3.1 
which  were  used  to  simulate  £>  (and  thus  the  p)  data 
may  be  meaned  to  form  l°x  1°  (and  2°x2°)mean  anomalies: 

(6.1)  \g'  -  E(Agj  cos$ffli)/ 1(003^) 


Where  the  AgJ  on  the  right  side  is  for  30'x30'  anomalies 
with  mean  latitude  0mi  and  Ag’  on  the  left  hand  side 
is  for  the  mean  anomaly  of  the  larger  l°x  1°  (and  2°x2°) 
blocks. 

The  1°  mean  residual  anomalies  are  shown  in  Figure 
6.1,  and  the  gradient  of  the  1°  anomaly  field  is  shown 
in  Figure  6.2  with  a  contour  interval  of  10  mgals.  We 
again  notice  that  the  anomaly  field  is  comparatively  smooth 
from  the  north  edge  up  to  12®  S. ,  becomes  rough  as  we  approach 
the  trench  which  starts  around  14®  s.  and  proceeds  south 
along  the  nominal  longitude  of  187?5  E, 

We  will  now  examine  the  prediction  of  two  1®  residual  anomalies 
#1  ,  -13°  >  <(>  >  -14®  ,  187®  <  X  <  188®  ,  +51  mgals  and 
#2  ,  -15®  >  < |>  >  -16®  ,  187®  <  X  <  188®  ,  -99  mgals.  The 
first  is  a  large  positive  anomaly  just  outside  the  trench  , 
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while  the  second  la  a  large  negative  anomaly  just  inside 
the  trench.  We  will  in  particular  test  the  effects  of 
the  extent  and  density  of  data  on  the  prediction  by  exam¬ 
ining  (S.l)  and  try  to  determine  a  mode  it  data  requirement 
for  economy  in  computer  time  consistent  with  satisfactory 
prediction.  Both  the  M/S  and  the  S/M  arcs  would  be  used 
to  confirm  the  results  ,  while  taking  into  account  the 
selection  of  data  points  as  shown  in  Figure  5.3b. 

Each  anomaly  was  predicted  separately  with  data  arranged 
symmetrically  for  that  anomaly.  This  has  been  done  throughout 
in  this  study.  Initially,  some  tests  were  made  to  examine 
the  correlations  between  adjacent  anomalies  by  predicting 
several  anomalies  together.  The  correlations  were  of  the 
order  of  0.2  to  0.6  as  already  reported  in  the  error  analysis 
by  Hummel  (1980,  e.g.  Table  4,  p.31).  Also,  as  the  prediction 
of  several  anomalies  together  by  a  common  data  set  did 
not  result  in  optimally  symmetric  arrangement  of  data  for 
the  individual  anomalies,  further  tests  for  correlations 
were  not  useful  in  this  simulation  study  where  we  could 
individually  examine  the  prediction  discrepancy  (5.1)  for 
each  anomaly. 


Table  6.1  shows  the  effect  of  the  extent  of  data  on 


the  prediction,  where  we  have  used  P  data  for  i|i  »0!5  , 

1°  ,  i?5  ,  2?5  ,  3?  5  .  The  predicted  anomaly  Ag*  ,  prediction 


discrepancy  Ae  ,  and  the  estimated  standard  deviation  of 
prediction  aAg'  are  tabulated  using  (2.9),  (5.1)  and 
(2.10).  The  results  are  given  only  for  Op  -  0.05  mgals, 
though  other  values  were  also  tried. 


Table  6.1  confirms  the  conclusions  from  Figure  5.6 
that  data  beyond  tj>  ■  1°  does  not  alter  the  predicted  value 
much;  in  fact  for  these  two  anomalies  the  prediction 
discrepancy  e  is  least  at  ■  1°  .  Secondly,  o^g*  is 
a  nominal  estimate  of  the  standard  deviation  of  prediction 
depending  primarily  on  the  number  of  data  points  and  their 
geometry.  a^_,  should  thus  be  interpreted  with  caution, 
e.g.  for  the  case  of  #2  anomaly.  The  prediction  discrepancy, 
e  ,  when  available,  is  a  better  criterion  for  judging  the 
accuracy  of  prediction.  Thirdly,- we  note  that  the  predicted 
value  of  #2  anomaly  is  very  discordant.  The  reasca  for 
this  will  become  clear  later  in  Section  7.  However,  Table 
6.1  does  show  that  in  this  case  also,  the  extent  of  data 
beyond  1°  does  not  alter  the  predicted  value  much. 

This  was  borne  out  by  other  tests  also  with  ij;  *4°,  5°,  6°. 


Table  6.2  shows  the  effect  of  the  density  of  data  points 
on  the  prediction.  The  results  are  given  for  ag  *  0.05 
mgals.  Apparently,  there  is  no  variation  in  the^prediction 
discrepancy  if  the  data  grid  is  more  dense  than  about  0?4x0?4. 
This  agrees  with  the  ideas  of  Siinkel  (1981,  p.91)  that 
there  is  no  improvement  in  prediction  error  estimation 
if  the  data  grid  density  is  increased  beyond  one-third 
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Table  6.1  Variation  in  the  Predicted  Mean  l°x  1°  Residual 
Anomaly  with  Variation  in  Extent  (\|f)  of  p  Da 

a..  =  0.05ragals  ,  o*  ■  lum/s 


Arcs 

No.  &  Spacing  of  pData 

l 

#  AX°  A<j>° 

#lAnom.  Ag'*+51  mgals 


#2  Anom.  Ag'=-99mgals 


Table  6.2  Variation  in  the  Predicted  Mean  l°x  1°  Residual 
Anomaly  with  Variation  in  Density  of  p  Data 

a-  a  0.05mgals  ,  a.  »  lum/s 


r 

Arcs 

No.  &  Spacing  of 

#  AX° 

p  Data 

A  d>° 

Predicted  Residual  Anom. 
„  (mgals) 

#2  Anom. 

Ag'=*-S 

>9  mgals 

1°.0 

N/S 

18  0.40 

0.41 

-34 

64 

12.8 

tt 

tt 

29 

0.28 

-33 

66 

12.4 

It 

S/N 

21 

0.41 

-38 

61 

12.4 

»t 

M 

30 

0.28 

-38 

61 

12.3 

0?5 

N/S 

16  0.20 

0.28 

-36 

63 

12.9 

ft 

ft 

31  it 

0.14 

-37 

62 

12.8 

S/N 

14  " 

0.28 

-36 

63 

13.0 

ft 

tt 

29 

0.14 

-37 

62 

12.8 

of  the  correlation  length  of  the  underlying  covariance  function, 
which  for  our  case  is  about  1?2  to  1?3  (Fig.  5.2  and  5.1 
respectively).  Our  results  are  also  in  general  agreement 
with  Rummel  et  al.  (1979,  Fig.  2,  p.  357),  who  found  the 
optimum  spacing  of  data  grid  at  250  km  altitude  to  be 
0?5x0?5  for  predicting  1°  anomalies. 

An  additional  test  was  tried  to  combine  the  extent 
and  density  of  data  by  having  a  0.2°  x  0.3°  grid  (AAxA<|>) 
till  ip  »  0?5  ,  a  0?4  x  0?4  grid  for  0? 5  <  \|><  1?0  ,  a  0?8x0?7 
grid  for  1?0<^<2?5  ,  but  it  did  not  lead  to  any  improvement 
over  a  grid  of  0?4  x  0?4  for  ^  “  1°  . 

We  thus  have  a  very  modest  data  requirement  leading 
to  economical  data  processing.  With  p  data  and  location 
of  satellites  on  tape,  several  anomalies  could  be  predicted 
one  at  a  time  in  the  same  run  with  data  being  selected 
optimally  for  each  anomaly  as  in  Figure  5.3b.  The  average 
run  time  on  Amdahl  470  V/6-II  computer  was  12  seconds  for 
predicting  each  l°x  1°  mean  anomaly.  This  time  due  to  the 
rigorous  computation  of  the  covariances  could  be  considerably 
reduced  by  interpolation  from  stored  tables,  as  discussed 
in  the  beginning  of  Section  5. 


7.  Prediction  of  l°x  1°  Mean  Residual  Anomalies 


The  optimum  data  configuration  from  Sections  5  and  6 
for  predicting  1°  mean  anomalies  is  then  to  use  p  data 
in  a  symmetric  grid  of  0?4  x  0?4  for  ip  »  1°  around  a  point 
which  is  down  the  arc  from  the  center  of  anomaly  block 
by  a  distance  equal  to  one-half  the  separation  between 
the  satellites. 

An  explanation  for  the  large  prediction  discrepancy 
for  #2  anomaly  in  the  trench  (see  Fig.  6.2)  may  be  the 
very  sharp  gradient  of  the  anomaly  field.  If  we  compare 
the  anomaly  field  in  Figure  6.1  with  the  accelerations 
in  Figure  5.5,  it  is  apparent  that  the  satellite  is  sensing 
the  "highs"  in  the  anomaly  field  of  51  ragals  at  13?5  S., 
187?5  E.  (# 1  anomaly)  and  that  of  77  mgals  at  17?5  S., 

186?5  E.  The  intervening  feature  of  H 2  anomaly  at  15?5  S., 
187? 5  E.  is  being  sensed  in  Figure  5.5,  but  is  considerably 
smoothed  out  because  of  the  two  "highs"  on  either  side. 

To  confirm  this  presumption,  we  may  examine  several 
anomalies  in  the  comparatively  smooth  field  away  but  near 
the  trench  and  then  extend  this  sample  to  the  rough  anomaly 
field  as  we  approach  the  trench  and  also  consider  anomaly 
blocks  over  the  trench.  Table  7.1  shows  24  predicted 
anomalies  in  the  area  bounded  by  latitudes  8°  to  12°  S. 
and  longitudes  184°  to  190°  E.  The  p  data  density  is 
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0?4x0?4  with  extent  if>  =>  1°  ,  and  the  results  are  for  o3 
as  0.05  ragals.  The  data  is  selected  symmetrically  individ¬ 
ually  for  each  anomaly  as  shown  in  Figure  5.3b. 


The  anomalies  are  listed  in  N/S  tiers  by  longitude 
to  follow  the  pattern  of  p  data  sensed  by  polar  satellites. 
The  RMS  A(root  mean  square)  value  of  anomalies  Ag'  ,  predicted 
value  Ag'  ,  and  the  prediction  discrepancy  e  is  also 
listed  in  sets  of  eight  anomalies.  The  minimum,  maximum 
and  mean  value  of  e  is  also  listed  for  sets  of  eight 
anomalies.  We  also  list  the  correlation  coefficient  r  A 
between  the  "true"  value  Ag'  and  the  predicted  value  Ag1 
of  the  1°  mean  anomalies: 


(7.1) 

where 


I  (Ag*  •  Ag' )/n 
i=l 


n  i  n  , 

(  l  Ag' 2/n)*(  V  Ag' 2/n)* 
i=l  i=l 

n  is  the  number  of  anomalies  being  compared. 


We  thus  find  very  good  agreement  in  this  area  between 
the  input  (1°  anomalies  obtained  from  input  30’  anomalies) 
and  the  predicted  1°  anomalies.  For  all  24  anomalies, 
the  RMS  value  for  the  input  Ag'  ,  predicted  Ag'  ,  and 
the  prediction  discrepancy  e  are  9.5,  6.2  and  5.5  mgals 
respectively  with  a  correlation  coefficient  r  of  0.84. 

The  nominal  estimate  of  the  standard  deviation  of  12.5 
mgals  is  therefore  pessimistic  by  a  factor  of  about  2.3. 

This  appears  to  be  due  to  the  use  of  a  global  covariance 
function  in  equation  (2.13)  instead  of  a  local  covariance 
function,  because  the  former  represents  onl”  the  average 
accuracy  of  the  estimated  quantities.  It  xs  well  known 
that  the  predicted  quantity  is  not  very  sensitive  to  changes 
in  the  covariance  function,  but  the  estimate  of  standard 
deviation  is  strongly  dependent  on  the  covariance  function. 
Schwarz  and  Lachapelle  (1980,  p.33)  found  that  the  standard 
errors  from  the  global  covariance  function  of  Tscherning- 
Rapp  (equation  (2.13))  were, about  three  times  larger  than 
they  should  be  for  Canada.  Schwarz  (1976,  p.13)  points 
out  that  for  local  applications,  it  is  necessary  that  the 
local  covariance  function  should  fit  the  degree  variances 
of  the  local  anomaly  gradients.  This  was  not  the  case 
for  the  covariance  .  function  used  in  this  study. 


We  now  extend  the  prediction  to  the  rough  anomaly 
field  as  we  approach  and  then  enter  the  trench  area.  From 
Figures  6.2  and  6.1,  we.  chose  eight  anomalies  bounded  by 
the  latitudes  12°  to  16°  S.  and  longitudes  186®  to  188® 

E.  where  the  eastern  anomalies  are  along  the  axis  and  then 
directly  over  the  trench,  and  the  western  anomalies  are 
immediately  by  the  side  of  the  trench.  The  particulars 
of  prediction  are  listed  in  Table  7.2. 
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Table  7.1  Prediction  of  l°x  1°  Mean  Residual  Anomalies  in 
the  Area  8aT  12°  S.,  184°-  190°  E. 

Extent  of  P  >/ata  V  -  1°  ; 

Density  of  Data  0?4x0?4  ; 

a..  »  o.05mgals  ,  a,  ■  lum/s 


Mean 


or  all  24  R 
anom.  Mean 


0.84 


Table  7.2  Prediction  of  l°x  1°  Mean  Residual  Anomalies  in 
the  Area  12°-16°S.  ,  186°  to  188°  E. 

Extent  of  P  Data  »  1°  ; 

Density  of  Data  0?4  x  0?4  ; 


o»  »  0.05mgals  ,  ■  lum/s 


Center  of 

Prediction  of  Anomaly 

Statist ics 

Correln . 

seq 

Anom. 

Block 

(mgals) 

for  e 

Coeff . 

ft 

d>° 

Yc 

Ag' 

Ag' 

e 

°Ag' 

(mgals) 

r 

1 

-12.5 

186 . 5 

-28 

-9 

19 

12.5 

2 

-13.5 

II 

-20 

-20 

-1 

II 

Min .-38 

3 

-14.5 

II 

-43 

-32 

11 

II 

Max.  19 

4 

-15.5 

II 

51 

13 

38 

II 

Mean  -2.2 

0.87 

MS 

37'.  3 

20.4 

21.9 

Mean 

-9.9 

-12.1 

-2.2 

12.5 

5 

-12.5 

1§7.5 

-<23 

6 

28 

i3.5  1 

6 

-13.5 

ti 

51 

57 

6 

12.4 

Min.  6 

7 

-14.5 

II 

-34 

18 

52 

12.5 

Max.  71 

8 

-15.5 

" 

-99 

-28 

71 

" 

Mean  39.3 

0.63 

5m§" 

59.2 

33.0 

46.3 

Mean 

-26.2 

13.1 

39.3 

12.5 

The  reasons  for  the  large  prediction  discrepancies  are 
clear  from  Figure  7.1,  where  the  "true"  anomalies  (implied 
by  30'  anomalies  used  for  simulating  p  data)  are  joined 
by  full  lines  and  the  predicted  anomalies  are  joined  by 
dashed  lines.  The  anomalies  are  shown  in  longitudinal 
profiles  from  184? 5  to  189?  5  E.  and  cover  the  area  from  8°S.  to 
16°  S.;  i.e.  all  anomalies  listed  in  Tables  7.1  and  7.2. 

Two  additional  anomalies  have  been  shown  on  each  side  of 
the  profiles  al  ng  the  longitudes  186? 5  and  187? 5  E.  We 
recall  that  the  profile  along  187?5  E.  goes  directly  over 
the  trench  south  of  latitude  14°  S.,  while  the  profile 
along  186? 5  E.  is  just  by  the  side  of  the  trench. 

We  see  from  the  above  two  profiles  that  the  satellite 
clearly  senses  the  "highs"  and  "lows"  of  the  anomaly  field, 
but  the  intervening  features  are  smoothed.  Notice  particularly 
the  anomalies  along  longitude  186?5  E.  with  centers  at 
latitudes  10?5  S. ,  14?5  S.  and  16?5  S.  which  have  small 
prediction  discrepancy  but  the  intervening  anomalies  at 
12? 5  S.  and  15? 5  S.  have  large  prediction  discrepancy  because 
the  satellite  does  not  sense  the  rapid  variation  in  the 
anomaly  gradient.  Note  that  it  is  not  the  magnitude  of 
the  anomaly  or  the  large  anomaly  gradient  which  causes 
prediction  discrepancy  (e.g.  the  anomalies  with  centers 
14?5  S.  and  16?5  S.  and  17?5  S.)  but  the  rapid  variation 
in  the  anomaly  gradient. 


Figure  7.1  Prediction  of  l°x  1°  Mean  Residual  Anomalies  (mgals).  Og«0.05mgals 
o.*lum/s.  "True"  Anomalies  (implied  by  30'  anomalies  used  for  simulating 
p  data)  joined  by  full  lines.  Predicted  anomalies  joined  by  dashed  lines. 
Extent  of  data  iji  ■  1°  ;  Density  of  data  0?4  x  0?4 


This  is  again  seen  in  the  profile  along  longitude  187?5  E., 
where  the  anomalies  with  centers  at  latitudes  ll!5  and 
13?5  S.  have  small  prediction  discrepancy  but  the  anomaly 
at  12?5  S.  has  ?  large  discrepancy  due  to  a  negative  anomaly 
gradient  changirg  to  a  positive  anomaly  gradient  at  12?5  S. 

The  variation  in  the  anomaly  gradient  is  so  rapid  in  the 
trench  area  after  15®  S.  that  the  predicted  anomaly  is  nearly 
zero  at  16?5  and  17!5  S.  The  rapid  variation  in  the  1°  anomaly 
gradient  is  seen  in  Figure  6.2.  It  is  also  Instructive 
to  see  the  variation  of  p  sensed  at  the  satellite  in 
this  area  in  Figure  5.5b. 

To  confirm  the  presumption  that  it  is  the  variation 
in  the  anomaly  gradient  that  causes  prediction  discrepancy, 
and  that  it  is  not  caused  by  lack  of  resolution  at  l°x  1° 
block  size,  we  will  examine  the  prediction  of  30'x  30'  and 
2°x  2°  anomalies  in  Sections  8  and  9  respectively. 


8.  Prediction  of  30* x  30'  Hear  Residual  Anomalies 


Table  8.1  lists  the  predicion  of  32  30*  anomalies,  in  sets 

of  8  anomalies,  in  the  area  bounded  by  latitudes  8°  to  12°  S. 

and  longitudes  186°  to  188°  E.  Analogously  to  the  data 

selection  criterion  for  predicting  1°  anomalies,  we  choose 

the  extent  of  data  for  predicting  30*  anomalies  as  i|;  «0?5  , 

with  density  of  p  data  grid  as  0?2x0?3  (AXxA<J>).  The 

data  is  selected  individually  for  each  anomaly,  symmetrically 

around  a  point  down  the  arc  from  the  center  of  anomaly 

block  at  a  distance  equal  to  one-half  the  separation  betwe*.  n 

the  satellites.  The  results  in  Table  8.1  are  for  a..  *»0.05i  ;als. 

p 

For  32  predicted  anomalies,  the  RMS  values  for  the 
input  Ag1  ,  predicted  A^g'  ,  and  the  prediction  discrepancy  c 
are  12.2,  10.5  and  7.3  respectively.  Again,  the  nominal 
estimate  of  the  standard  deviation  of  18.1  mgals  is  pessi¬ 
mistic  by  a  factor  of  about  2.5.  The  input  and  the  predicted 
anomalies  have  been  plotted  in  Figure  8.1  along  four  longi¬ 
tudinal  profiles.  The  input  anomalies  have  been  joined 
by  full  lines  and  the  predicted  anomalies  by  dashed  lines. 

It  is  clear  that  the  trend  of  input  30'  anomalies  is  quite 
faithfully  followed  by  the  predicted  anomalies.  In  this 
area  of  comparatively  smooth  anomaly  field  (also  see  Fig. 

3.2),  we  generally  have  a  resolution  of  30*  anomalies  but 
it  is  perhaps  more  indicative  of  satisfactory  modeling 
of  the  "process"  and  the  "observational"  noise. 

The  average  run  time  for  the  prediction  of  a  30'  x30’ 
mean  anomaly  with  the  specified  data  grid  selection  was 
only  7  seconds  on  Amdahl  470  computer.  This  time  could 
be  reduced  further  if  the  covariances  were  approximated 
by  interpolation  from  stored  tables,  instead  of  being  computed 
rigorously . 
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Table  8.1  Prediction  of  Mean  30 'x  30'  Residual  Anomalies  in 
Area  8°  to  12°  S.,  186°  to  188°  E. 

Extent  of  p  Data  ^»0?5  ; 

Density  of  Data  0?2x0?3  ; 


O.OSmgals  , 


lum/s 


Center  of 
Anom.  Block 


-8.25 

-8.75 

-9.25 

-9.75 

•10.25 

-10.75 

•11.25 


-8.75 

-9.25 

-9.75 

•10.25 

•10.75 

•11.25 

-11.75 


186.25 


Prediction  of  Anomaly  (mgals) 

Ag'  *V  *  °Ag. 

Statistics 
for  e 

(mgals) 

-15 

-10 

5 

18. C 

-19 

-7 

13 

18.1 

Min.  2 

-6 

3 

9 

II 

7 

9 

2 

18.0 

Max. 13 

5 

11 

6 

18.1 

-2 

4 

6 

II 

Mean  6.7 

-10 

-1 

9 

It 

-6 

-1 

5 

" 

Coef f . 


-8.25 

186.75 

-8.75 

II 

-9.25 

II 

-9.75 

II 
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II 

-10.75 

11 

-11.25 

II 
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II 

RMS  10.2 
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-8 

25 
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-8 
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II 

-9 
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Figure  8.1  Prediction  of  30'x30'Mean  Residual  Anomalies  (mgals) 
o«  ■  0.05  mgals;  a»,*lum/s  Input  anomalies  used  for  simulating  p 
data  joined  by  full  lines.  Predicted  anomalies  .joined  by  dashed 
lines.  Extent  of  data  Density  of  data  0?2x0?3 


9.  Prediction  of  2*x  2°  Mean  Residual  Anomalies 


To  further  confirm  that  large  prediction  discrepancies  for 
1°  anomalies  near,  and  in  the  trench  area  are  due  to  the 
rapid  variation  in  the  anomaly  gradient  and  not  due  to  the 
resolution  limit  of  the  block  size,  we  examine  the  prediction 
of  2°x  2°  mean  anomalies.  The  2°  residual  anomalies  implied 
by  the  30 1  residual  anomalies  (equation  (6.1))  are  shown 
in  Figure  9.1,  and  the  gradient  of  this  anomaly  field  is 
shown  in  Figure  9.2  at  a  contour  interval  of  5  mgals.  If 
we  compare  Figure  9.2  with  the  $  accelerations  in  Figure 
5.5,  it  is  immediately  apparent  that  much  greater  resolution 
is  possible  from  the  accelerations  than  can  be  depicted 
by  2°  anomaly  field  (cf  Fig.  6.2  and  3.2  for  1°  and  30' 
anomaly  fields). 

Secondly,  we  notice  from  Figure  9.2  that  the  outline  of 
trench  from  14°  S.  proceeding  south  along  187? 5  E.  has  been 
considerably  generalized  due  to  the  coarseness  of  the  resolution 
of  2°  blocks  for  depicting  this  short  wavelength  feature. 

The  variation  in  anomaly  gradient  has  become  the  more  impor¬ 
tant  feature  affecting  three  2°  blocks  centered  at  17°  S., 

185°  E.;  17°  S.>  187°  E.  and  15°  S. ,  187°  E.  These  blocks 

would  have  large  prediction  discrepancy  because  of  the  rapid 
variation  in  anomaly  gradient. 

We  in  fact  find  this  to  be  so.  Figure  9.3  shows  the 
plot  of  "true"  and  predicted  anomalies  and  Table  9.2  lists 
the  prediction  of  15  2°  anomalies  in  the  area  bounded  by 
latitudes  8°  to  18°  S.  and  longitudes  184°  to  190°  E.  The 
extent  of  data  was  ip  ■  1®  for  each  individual  anomaly  sym¬ 
metrically  around  a  point  r  s  in  Figure  5.3b,  with  data  grid 
0?4°x0?4  ,  and  the  result?,  are  quoted  for  o«  »  0.05  mgals. 

If  we  omit  all  2°  blocks  south  of  latitude  14°  S.  from 
consideration  due  to  large^  variations  in  anomaly  gradients, 
the  RMS  value  of  "true"  Ag'  ,  predicted  Ag'  ,  and  prediction 
discrepancy  for  9  anomalies  is  6,1,  6.6,  5.3  respectively. 
Perhaps  a  larger  value  of  ip  ,  say  2°,  would  have  resulted 
in  better  prediction. 
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*M«n,ro  o  3  Prediction  of  2°x  2°  Mean  Residual  Anomalies  (mgals) 
Figure  .  OA-lum/s.  "True"  anomalies  (implieu  by  30 

to  »  i°  •  Density  of  data  O.Axo.e 
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Summary  and  Conclusions 


zX  The  purpose  of  this  simulation  study  was  to  test  the 
accuracy  of  prediction  of  mean  gravity  anomalies  from  the 
GRAVSAT  mission,  to  establish  the  resolution  achievable 
in  terms  of  anomaly  block  size,  and  to  draw  conclusions 
for  the  data  requirement  for  economical  processing  in  terms 
of  computer  time... 


Two  low  satellies  were  considered  in  identical  nominally 
polar  orbits  in  a  purely  gravitational  field  at  an  altitude 
of  150  km  over  a  spherical  earth,  and  with  an  along-track 
separation  of  200  km  between  the  satellites.  The  input 
gravitational  field  was  described  by  1500  (residual)  30'x30’ 
altimeter  derived  anomalies  (R.H.  Rapp,  private  communication, 
1980),  from  which  the  anomaly  implied  by  the  gravitational 
model  GEM9  had  been  subtracted  out.  The  anomalies  covered 
the  area  0°  to  25°  S.  in  latitude,  and  180°  to  195°  E. 
in  longitude  in  the  south-west  Pacific  Ocean,  including 
the  sharply  varying  anomaly  field  associated  with  the  Tonga 
Kermedec  trench,  which  lay  in  the  center  of  the  area  proceeding 
in  a  generally  southerly  direction  from  about  14°  S. 


23,128  range-rate  (p)  measurements  were  simulated  by 
Wong  and  Sjogren  (1980)  with  data  interval  of  2  seconds, 
(i.e.  A$  a  0?14)  and  the  separation  of  north  to  south,  as 
well  as  south  to  north,  arcs  each  being  0?2  in  longitude 
(i.e.  AX  *  0?2).  This  corresponds  to  a  mission  duration 
of  109.4  days  (Jekeli  and  Rapp,  1980,  p.12),  or  a  little 
less  than  4  months. 


The  p  measurements  were  impressed  with  a  random 
"observational"  noise  of  mean  0  and  standard  deviation 
lym/s  to  simulate  the  precision  of  the  GRAVSAT  mission, 
and  supplied  to  us  for  this  study.  The  "observations"  cor¬ 
respond  to  residual  range-rate,  i.e.  range-rate  observed 
in  the  earth's  gravitational  field  minus  the  computed  range- 
rate  in  GEM9. 


The  "noisy"  residual  range-rate  p  was  filtered  of 
the  observational  noise  by  approximating  J>  in  the  least 
squares  sense  by  a  cubic  spline  function  $ s  with  fixed 
knots.  It  was  determined  from  numerical  tests  that  a  spline 
knot  spacing  of  10  seconds  filtered  out  a  noise,  whose  RMS 
value  varied  from  0.92  to  1.04ym/s  for  various  arcs.  The 
first  derivative  of  the  "smoothed"  j5s  was  obtained  analytically 
from  the  spline  function,  providing  us  the  residual  line 
of  sight  acceleration,  p  (see  equations  2.4  and  2.5). 

p  values  were  used  in  this  study  as  the  "observations", 
as  these  could  be  easily  related  to  the  anomalous  potential 
at  the  satellite  locations.  Further,  the  procedure  of  obtaining 
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p  from  noisy  p  measurements  also  largely  filters  out 
any  systematic  errors  causing  long  wavelength  perturbations 
e.g.,  due  to  uncertainties  in  the  determination  of  initial 
state  vectors  of  the  satellites  (see  page  21). 

The  residual  anomalies  A* g'  were  predicted  in  block  sizes 
of  30'x30'  ,  l°x  1°  ,  2°x  2°  using  least  squares  collocation, 
as  in  equation  2.9: 

Ag'  “  A  Tjj,  ;  A  “  —  »  £  “  — T4Ti+  — 

where  Tjj,  =  g_  is  the  "observational"  data  used  in  the  pre¬ 
diction,  with  the  associated  noise  matrix  I)  .  is 

the  auto-covariance  matrix  of  the  p  data, “while  C^g'Tjj, 
is  its  cross-covariance  vector  with  the  residual  anomaly 
of  specified  block  size.  A  may  be  called  the  estimating 
operator  or  weight  function  (see  Fig.  5.6,  p.26),  which 
may  be  examined  to  judge  if  distant  data  may  be  excluded 
without  significatly  affecting  the  predicted  value.  The 
underlying  covariance  functions  Ctjits,  and  cAg'T£  may 
also  be  examined  for  their  correlation  distanced  see  Fig. 

5.1,  5.2,  pp.  17,  19). 

It  was  inferred  from  these  considerations  that  for 
predicting  residual  l°x  1°  mean  anomalies  (residual  to  GEM9, 
complete  to  degree  and  order  20),  a  data  cap  of  spherical 
radius  1°  is  adequate.  More-  distant  data  not  only  contribute 
to  the  prediction  with  decreasing  amplitude,  but  also  with 
varying  sign  tending  to  largely  cancel  out  the  contribution. 
Further,  there  is  no  improvement  in  the  predicted  , value 
of  1°  residual  anomalies  if  density  of  data  is  increased 
beyond  0?4x0?4  ,  which  is  about  one-third  of  the  correlation 
distance  of  .  These  conclusions  were  borne  out 

by  numerical  tests  in  Tables  6.1  and  6.2  by  examining  the 
prediction  discrepancy  e  in  equation  5.1: 

e  =  tfg'  -  Ag ' 


where  Ag'  is  the  predicted  value  and  Ag'  is  the  input, 
or  the  implied  residual  anomaly  obtained  by  meaning  the 
input  30'  residual  anomalies  used  for  simulating  5  measure¬ 
ments. 

This  leads  to  a  very  modest  data  requirement  for  predicting 
1°  anomalies,  i.e.  a  data  density  of  about  0?4x0?4  in  a 
data  cap  of  spherical  radius  1°  ,  a  total  of  about  20  data 
points.  The  needed  density  would  be  available  in  a  GRAVSAT 
mission  of  about  one  month.  Longer  duration  of  mission 
would  permit  detection  of  systematic  errors  and  editing 
of  data  of  doubtful  quality 
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Several  1°  anomalies  may  be  predicted  in  a  single  run, 
separately  selecting  data  symmetrically  for  each  individual 
anomaly.  The  average  time  for  predicting  each  1°  anomaly 
was  about  12  seconds  on  Amdahl  470  V/6-II  computer.  In  this 
study,  the  covariances  were  rigorously  computed  necessitating 
16  calls  to  COVAX  (Tscherning,  1976)  for  computing  each 
element  of  CtAt*  *  and  similar  lengthy  computations  for 
C/c'Tj  •  The  time  for  prediction  of  each  1°  anomaly  could 
be  substantially  reduced,  by  interpolating  covariances  from 
stored  tables. 

It  was  noted  from  Fig.  5.2  and  5.6  for  and 

the  estimating  operator  A  that  the  data  cap  Should  be 
centered  not  at  the  center  of  anomaly  block  being  predicted, 
but  at  a  distance  from  the  block  center  down  the  arc  equal  to 
one-half  the  separation  between  the  satellites.  This  is 
also  shown  in  Fig.  5.3b  (page  20). 

Rummel  (1980,  p.13)  had  related  ojj  ,  the  noise  of  p 
as  0 .05  mgals  corresponding  to  the  noise  of  £  of  lym/s 
at  satellite  altitude  of  150  km.  This  was  numerically  confirmed 
in  this  study  by  examining  the  estimating  operator  A  (p.26) 
for  different  values  of  ox  .  For  smaller  values  than  0.05 
mgals,  the  oscillations  of  A  were  with  much  greater  amplitude 
resulting  from  unstable  inversion  of  CtaTJL  •  On  *^e  0l-^er 
hand,  for  values  of  o^  larger  than  0.0b  mgals,  the  oscil¬ 
lations  of  A  were  damped  out. 

24  l°x  1°  residual  anomalies  were  predicted  in  the  area 
8°  to  12°  S.  latitude  and  184°  to  190°  E.  longitude,  which 
was  near  but  away  from  the  trench,  and  thus  with  a  comparatively 
smooth  anomaly  gradient  field  (p.29).  The  RMS  prediction 
discrepancy  for  the  24  anomalies  was  5.5  mgals,  while  the 
RMS  value  of  "true"  anomalies  was  9.5  mgals.  The  correlation 
coefficient  between  the  predicted  and  true  anomalies  was 
0.84  (see  page  34) . 

This  showed  a  good  agreement  for  "smooth"  anomaly  gradient 
field.  The  prediction  of  1°  anomalies  was  then  tried  in 
the  "rough"  anomaly  gradient  field  while  approaching,  and 
then  directly  over  the  trench.  This  was  done  in  latitude 
12°-18°  S.  along  two  longitudinal  profiles  of  186®5  and 
187?5  E.  The  first  profile  was  immediately  to  the  west 
of  the  trench,  while  the  second  profile  was  directly  over 
the  trench.  It  was  found  (page  36)  that  though  the  'highs 
and  "lows"  in  the  anomaly  field  were  generally  predicted 
satisfactorily,  the  intervening  anomalies  were  smoothed 
out  between  the  extreme  values.  Further,  the  prediction 
was  poor  whenever  the  sign  of  anomaly  gradient  changed; 
and  particularly  over  the  trench,  the  very  rapid  short  wavelength 
variations  in  the  anomaly  gradient  were  practically  sensed 
as  zero  at  the  satellite  altitude.  This  is  also  seen  in 
the  plot  of  accelerations,  p  ,  sensed  at  the  satellite 
(see  page  25). 


The  question  of  the  resolution  achievable  in  terms  of 
anomaly  block  size  was  addressed  by  predicting  32  residual 
30'x30'  anomalies  in  the  smooth  anomaly  gradient  field  in 
the  area  8°  to  12°  S.  latitude  and  186°  to  188°  E.  longitude. 

The  RMS  prediction  discrepancy  for  the  32  anomalies  was 
7.3  mgals,  while  the  RMS  value  of  "true"  anomalies  was  12.2 
mgals  (see  page  38).  The  plot  of  predicted  and  true  (i.e. 
input)  30*  anomalies  shows  (see  page  39)  that  while  the 
general  trend  of  30*  anomalies  was  faithfully  predicted, 
the  individual  short  wavelength  feature  in  30'  anomalies 
were  not  correctly  predicted  particularly  when  the  sign 
of  anomaly  gradient  changed  sharply. 

The  problem  of  prediction  when  anomly  gradients  are 
changing  sharply  was  again  found  while  predicting  2°x  2° 
anomalies  in  the  area  8°  to  18°  S.  in  latitude,  and  184° 
to  190°  E.  in  longitude  (page  44).  9  residual  2°  anomalies 
between  8°  to  14°  S.  latitude  had  a  RMS  prediction  discrepancy 
of  5.3  mgals,  while  the  RMS  value  of  true  anomalies  was 
6.1  mgals.  The  prediction  in  the  sharply  varying  anomaly 
gradient  field  south  of  14°  S.  was  unsatisfactory  (page  42). 

It  was  however  clear  that  there  was  considerably  more  detail 
in  the  accelerations  sensed  at  the  satellite  (page  25)  compared 
to  the  variations  which  may  be  represented  in  the  2°  anomaly 
field  (page  42). 

It  was  also  found  that  the  formal  estimate  of  the  standard 
deviation  of  the  prediction  (equation  2.10)  was  2.3  to  2.5 
times  larger  then  the  RMS  prediction  discrepancy.  This 
was  due  to  the  use  of  a  global  covariance  function  (Tscherning 
and  Rapp,  1974)  instead  of  a  local  covariance  function, 
which  may  have  permitted  a  better  fit  to  the  variance  of 
anomaly  gradients  in  the  area.  This  parallels  the  conclusion 
of  Schwarz  and  Lachapelle  (1980),  who  found  the  estimate 
of  standard  deviation  of  the  prediction  from  the  global 
covariance  function  to  be  about  three  times  larger  in  Canada. 

We  may  conclude  that  this  simulation  study  has  established 
satisfactory  procedures  for  predicting  the  anomalous  gravi¬ 
tational  field  from  the  range-rate  data  from  GRAVSAT  mission. 

It  has  also  established  the  requirement  of  data  for  economical 
processing  in  terms  of  computer  time.  The  resolution  achievable 
is  l°x  1°  anomaly  blocks  with  an  accuracy  of  prediction  of 
about  5.5  mgals.  But  this  is  possible  only  in  areas  where 
the  anomaly  gradients  are  not  varying  sharply,  e.g.  areas 
associated  with  trenches.  Near  the  trench,  the  anomaly 
highs  and  lows  may  be  predicted  satisfactorily,  but  the 
intervening  anomalies  would  be  smoothed  out  between  the 
extreme  values.  However,  directly  over  the  trench,  the 
prediction  from  range  rate  data  at  satellite  altitude  would 
be  unsatisfactory  due  to  the  very  rapid  short-wavelength 
variations  in  the  anomaly  gradients,  which  is  not  resolved 
at  the  satellite  altitude. 
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~7  It  was  found  fa  this  studyv)that  though  the  estimates 
of  the  formal  standard  deviation  of  the  predicted  1*  anomalies 
were  reduced  with  the  Increase  in  size  of  the  data  cap, 
and  with  the  increase  in  the  density  of  data  inside  the 
cap,  there  was  1J ctle  improvement  in  the  predicted  value 
of  1 6  anomalies  by  least  squares  collocation  technique. 

Besides,  there  was  a  large  increase  in  the  computer  time 
required  for  the  processing  of  the  increased  data--  Further 
studies  are  required  to  determine  if  other  techniques,  or 
a  modification  of  least  squares  collocation  technique,  would 
allow  the  economical  processing  of  increased  data.  And\ 
if  the  use  of  the  increased  data  would,  in  addition  to  lowering 
the  estimate  of  the  standard  deviation,  \lso  result  in  signifi¬ 
cant  improvement  of  the  predicted  value  of  the  anomalies. 
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